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ABSTRACT 

This paper presents the detection, interpretation and removal of the signal resulting from interactions of high energy particles with the Planck 
High Frequency Instrument (HFI). These interactions fall into two categories, heating the 0.1 K bolometer plate and glitches in each detector 
time stream. Glitch shapes are not simple single pole exponential decays and fall into three families. The glitch shape for each family has been 
characterized empirically in flight data and removed from the detector time streams. The spectrum of the count rate per unit energy is computed for 
each family and a correspondence to where on the detector the particle hit is made. Most of the detected glitches are from galactic protons incident 
on the Si die frame supporting the micromachinced bolometric detectors. At HFI, the particle flux is ~ 5 cirT 2 s and is dominated by protons 
incident on the spacecraft with an energy >39 MeV, leading to a rate of typically one event per second and per detector. Different categories of 
glitches have different signature in timestreams. Two of the glitch types have a low amplitude component that decays over nearly 1 second. This 
component produces an excess noise if not properly removed from the time ordered data. We have used a glitch detection and subtraction method 
based on the joint fit of population templates. The application of this novel glitch removal method removes excess noise from glitches. Using 
realistic simulations, we find this method does not introduce signal bias. 

Key words. Cosmology: observations - Cosmic background radiation - Surveys - Space vehicles: instruments - Instrumentation: detectors - 
Space radiatiomparticle environment 



1. Introduction 

This paper, one of a set associated with the 2013 release of data 
from the Planck 1 mission (Planck Collaboration I 2013), de- 
scribes the detection of high energy particles with the 0.1K in- 

1 Planck (http://vjww.esa.int/Planck) is a project of the 
European Space Agency (ESA) with instruments provided by two sci- 
entific consortia funded by ESA member states (in particular the lead 
countries France and Italy), with contributions from NASA (USA) and 



strument and bolometric detectors in the Planck High Frequency 
Instrument (HFI)(Lamarre et al. 2010) and removal of these sys- 
tematic effects from the millimeter and sub-millimeter signals. 
Bolometers (Holmes et al. 2008), such as those used in HFI 
on Planck, are phonon mediated thermal detectors with finite 
response time to changes in the absorbed optical power. The 



telescope reflectors provided by a collaboration between ESA and a sci- 
entific consortium led and funded by Denmark. 
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Fig. 1. Top left and right: Completed multimode 545GHz and 
single mode 143GHz SWB bolometer modules. Middle: An ex- 
ploded view of the assembly of a PSB. Alignment pins, shown in 
solid black, fix the aft and fore bolometer assemblies to an angu- 
lar precision of < 0.1°. The SWB assembly is similar to the PSB 
aft bolometer assembly and does not include a feedhorn aperture 
integrated with the module package. Bottom: A PSB pair epox- 
ied in the module parts prior to mating. To the right, the feedhorn 
aperture can be seen through the fore bolometer in the housing. 
To the left, the 1/4 wave backshort can be seen through the aft 
bolometer absorber mesh. 



bolometer consists of a micromachined silicon nitride (Si-N) 
mesh absorber with a Ge thermistor suspended from the Si die 
frame. Each bolometer is mounted to a metal housing which 
are assembled into two types of bolometer modules shown in 
Figure 1; a Spider Web Bolometer (SWB) which detects total 
power and a Polarized Sensitive Bolometer (PSB) (Jones et al. 
2007). In the PSB, there are two bolometers to measure indepen- 
dently power in each of the two linear polarizations, as shown in 
Figure 1 . The bolometers are mounted to a copper-plated stain- 
less steel plate cooled to 0.1K and stabilized to a few microK 
of the temperature set point. Two 'dark' bolometer modules are 
blanked off at the 0. IK plate and are used to monitor systematic 
effects. The 0.1 K bolometer plate is surrounded by a ~ 1. 5mm 
thick aluminum box cooled to 4.5 K. Light is coupled into each 
bolometer module using a feedhorn and filter stack at 4.5 K and 
1.5 K respectively (Lamarre et al. 2010). 

At L2, high energy particles from the sun and a galactic 
sources, primarily, protons, electrons and helium nuclei, are in- 



cident on the spacecraft. This particle flux causes two main ef- 
fects in HFI that have been reported previously. There is a time 
variable thermal load on the 0. IK plate (Planck Collaboration II 
2011) and a significant rate of glitches in the bolometric signal 
(Planck HFI Core Team 2011). In this paper, we report on the 
evolution of these effects over the entire mission, from July 3, 
2009 through January 14, 2012 and on the analysis technique 
used to remove the effect of glitches from the data. Three fam- 
ilies of glitches, long, short and slow (or longer as named in 
Planck HFI Core Team 2011), were found by comparing and 
stacking many events. We present the coincident counts, energy 
distribution, total counts and variations of glitch shapes within 
each family. These data taken together allow identification of 
the physical cause of the glitch events for two of the types. Long 
glitches are energy absorption events in the silicon die and short 
glitches are events in the optical absorbing grid or Ge thermistor. 
Experimental evidence obtained using beam line and lab tests us- 
ing spare HFI bolometers supports these conclusions (Catalano 
et al. 2013; Planck Collaboration II 201 1). The cause of the slow 
glitches, however, has not been identified. The new analysis pre- 
sented here employs a joint fitting and removal of these three 
glitch types. We show that this technique improves the noise 
performance of HFI. We report on rare effects not previously 
reported, including the response of the instrument to solar flares, 
secondary showers and large high energy events. We have sim- 
ulated the effect of the glitches on the data quality and find that 
systematic biases are small, < 10~ 4 compared to the cosmologi- 
cal signal. 

2. Glitch characterisation from bolometer data 

As reported previously (Planck HFI Core Team 201 1), the glitch 
rate in Planck HFI is about 1 per second. To remove glitches 
from the data time streams of many instruments with bolomet- 
ric detectors, (Crill et al. 2003; Masi et al. 2006), the instrument 
transfer function is deconvolved from the data. This has the ben- 
efit of increasing the signal to noise on each glitch which mini- 
mizes the fraction of flagged date. The populations of glitches in 
Planck HFI is unusual. Three distinct families of glitch transfer 
functions, shown in Figure 2, are found in the data. The three 
families of events, "short", "long" and "slow" are differentiated 
by their decaying time constants and amplitudes. Of these, the 
long glitches dominate the total rate. This type, however, has a 
significant amplitude over more than 1 second. Therefore, sim- 
ple flagging and removal of these glitches would lead to an un- 
acceptably large loss of data. To avoid such data loss, precise 
glitch templates, need to be subtracted from the time ordered 
data. To do this, we take advantage of the excellent pointing 
accuracy (Planck Collaboration I 2011) and redundancy of the 
scan strategy and have adopted a method that iterates between 
estimates of the signal from the sky and signal from the glitches 
(Planck HFI Core Team 201 1). A key to the effectiveness of the 
method is a comprehensive understanding of the three families 
of glitches and stability of these families over the course of the 
mission. In this section, we describe the general features of tem- 
plates of each of the three glitch families obtained from stacking 
many normalized glitches and the general features of the glitch 
templates derived from fitting the normalized decay profile with 
the sum of 4 exponentials Planck Collaboration VII (2013). 

2.1. Short glitches 

The short events have a rising time constant well within the sam- 
pling period of ~ 5 ms and then a fast decay. The highest am- 
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Fig. 2. Different types of glitches found in bolometer data. 
Examples for three distinct families of glitch transfer functions 
displayed for a typical the PSB-a bolometer. Events typical of 
the blue curves are called short events, those typical of the black 
curves are called long, and those of the red are called slow. 
Typical variation is in the shape within each family is shown 
for short and long glitches. 



plitude glitches are short glitches. The short glitch template for 
each bolometers is shown in Figure 3. We see small variability 
in the averaged short glitch template depending on the detec- 
tor. The shortest exponential decay has a time constant similar 
to the fast part of the optical transfer function with the three 
additional terms, or tail, having a very low amplitude < 1CT 2 
relative to the peak and intermediate time constants of typically 
40ms, 400ms and 2 seconds. The 2-second decaying signal is 
~ 10~ 4 and was detected only after stacking short glitches mea- 
sured over the entire mission. We have observed some similar- 
ities between the transfer function of short glitches and the op- 
tical transfer functions (Catalano et al. 2013). The amplitude of 
the tail is variable from event to event, but the amplitude of the 
tail is small enough such that its variability does not affect the 
processing of data. Moreover, we have found the existence of 
essentially two sub-categories of short glitches, one of each is 
shown in Figure 2, that can be distinguished by the amplitude 
of the tails. In fact, the tail is smaller for higher energy events 
and larger for lower energy events with a variation of a factor 
2. The time constants, as well as amplitudes resulting from the 
fit of three exponentials on short glitch measurements are given 
in Figure 4 for all bolometers. The fit has been performed after 
stacking events with amplitudes between 1200 and 2400 sigmas 
of the noise, with corresponds energies of the transition between 
the two sub-categories of short events. The 2 second time con- 
stant is not included in the template fit. We observe some scatter 
of the values across bolometers which are reflecting more the 
degeneracies in the fit of the different time constants and ampli- 
tudes than variations of the template shape between bolometers. 
A careful study of degeneracies of parameter fitting is performed 
in Planck Collaboration VII (20 1 3) for the transfer function anal- 
ysis. 



2.2. Long glitches 

The long glitches are on average lower amplitude than short 
glitches but dominate the rate over other categories of glitches 
by more than an order of magnitude. The rising time constant 
is very fast and within the sampling period as for short glitches. 
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Fig. 3. Short glitch templates for all detectors. Blue color is for 
PSB-a, green for PSB-b and violet for SWB. 
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Fig. 4. Parameters of the glitch templates for all bolometers. 
Left: short glitch template, right : long glitch template. The type 
of bolometers is indicated by color (blue = PSB-a, green=PSB- 
b, magenta= SWB) and the type of glitch by symbol (circle = 
long glitches, stars = short ones). Values quoted are obtained af- 
ter fitting exponentials on templates estimated by stacking large 
number of events. A three-point filter was applied to the data 
prior to the fit of exponentials. There are degeneracies in the 
fit of the different parameters, this is studied in detail in Planck 
Collaboration VII (2013) for the similar problem of transfer 
function fitting. 



The shortest exponential decay has a time constant similar to the 
fast part of the optical transfer function, followed by a tail with 
a much larger amplitude than for the short glitches. In Figure 5 
long glitch templates are shown for several bolometers. As with 
the short glitches they are estimated by the stacking of large 
number of events and fit to a sum of 4 exponential terms. 

The slow tail has typical time constants of 50 ms, 500 ms and 
2 seconds for PSB-a bolometers and 35 ms, 500 ms and 2s for 
PSB-b and SWB bolometers. The intermediate time constant of 
35 to 50 ms starts with an amplitude relative to the peak of «6% 
for PSB-a and of *3% for PSB-b and SWB bolometers. The 
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Fig. 5. Long glitch templates for all detectors. Blue color is for 
PSB-a, green for PSB-b and violet for SWB. 

long time constants of «500 ms and 2 seconds have an amplitude 
which is about 0.2% of the event peak. Given that most of the 
glitches are these long events in HFI data, the tail would cause 
excess power in the noise if not carefully removed. 

The most striking feature in in Figure 5 is that the intermedi- 
ate time constant is significantly slower and lower in amplitude 
for PSB-b and SWB as compared to PSB-a. The component of 
the tail with the longest time constant, of the order of 2 s, has 
a similar amplitude for all bolometers. The time constants and 
amplitudes are given in Figure 4 for all bolometers. 

We do not observe significant change in the time constants 
of the long glitch tail from glitch to glitch. However, the relative 
amplitude of the long tail with respect to the fast part show very 
significant dispersion. This is illustrated in Figure 6, where we 
show the distribution of the peak glitch amplitude as a function 
of the relative amplitude of the glitch tail for a typical PSB-a 
and PSB-b. If the long glitch shape was identical for all glitches, 
we would observe a distribution along a vertical axis centered at 
unity, with a width fixed by the measurement errors. Instead, we 
observe that the width of the distribution is larger than expected 
from just increasing signal to noise for increasing glitch peak 
amplitude. This intrinsic scatter in the long tail amplitude can 
only be partly due to the signal averaging on the sample affecting 
mainly the peak of the glitch, and decreasing the amplitude by a 
factor which depend on the timing of the glitch impact with in a 
sample period (Lamarre et al. 2010). Moreover, the distribution 
of slow part amplitude is for some bolometer double-peaked. 
Also, a small proportion of events have an amplitude of the fast 
part relative to the slow part which is different than for most 
long glitches. In particular, slightly less than 10% of events have 
a smaller amplitude tail than long glitches and about 0.5% have 
a higher tail. The events close to zero in the x-axis, are a small 
quantity events with low enough amplitude that short and long 
event cannot be distinguished. 

Moreover, we observe that the ratio between the amplitude 
of the slow tail and the peak amplitude of long glitches, is not a 
constant with respect to the peak amplitude of glitches. Higher 
energy events have a relatively higher longer tail. This non- 
linearity is seen in Figure 6, as a slight curvature of a the subset 
of high energy events highlighted by a dashed line. We have cor- 
rected for this nonlinearity by modelling the peak amplitude as a 
quadratic function of the tail amplitude. By comparing the non- 
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Fig. 7. Comparison of a long glitch (black) and a slow glitch 
(dot dashed blue). For this plot, the two high energy events were 
rescaled to match after 200 ms. The difference is shown in blue. 



linearity functions between bolometers, we observe a system- 
atic difference between PSB-a and PSB-b/SWB bolometers. The 
non-linearity for PSB-a corresponds to a factor 2 between the 
lowest and the highest energies of detected events, and is slightly 
stronger than that PSB-b and SWB which is a factor 1.5. Apart 
from this non-linearity, which is accounted in the glitch subtrac- 
tion method, we have checked that there is no dependency of the 
time constants of the slow tail with respect to the energy of the 
events, as illustrated by Figure 25 showing the stacking of events 
for different energy windows. 

2.3. Slow glitches 

For PSB-a bolometers only, the distribution of glitch peak ampli- 
tude to tail amplitude has a significant extra peak at ai/g max ~ 6. 
This extra peak is the signature of the slow glitches. The slow 
glitches have a rising time comparable to the optical time con- 
stant; much slower than the 2 other glitch families. The tail of the 
slow glitch is similar to that for a long glitch. In Figure 7 shows 
the comparison of a high amplitude long glitch event with a slow 
glitch. The intermediate time constant which is of the order of 50 
ms for long glitches is slightly but significantly larger for slow 
glitches as can be seen Figure 7. The two tails are proportional 
to a very good accuracy after m 200 ms after the peak ampli- 
tude. Even with these differences, the long glitch template for a 
given bolometer is a good proxy of the slow glitches in the same 
bolometer. 

2.4. Evolution of the Glitch Rate over the Mission 

Over the course of the whole mission, the rate for each type of 
glitch decreases, as shown in Figure 8 This decrease is universal 
for all bolometers as shown in Figure 9. The signal from diode 
sensor TC2, the most shielded of the 3 diodes in the Standard 
Radiation Event Monitor (SREM, Mohammadzadeh et al. 2003) 
located on the sun side of the PZancfcspacecraft, follows a simi- 
lar trend with time as shown in Figure 10. This trend is predicted 
for galactic cosmic rays modulated by the heliosphere of the sun 
(Gleeson & Axford 1968; Bobik et al. 2012) and is observed by 
ground stations on earth and other spacecraft in the solar system. 
(Usokin et al. 2011; Wiedenbeck et al. 2005; Adriana 2011). 
Taken together, the glitch rate, SREM diode signal and data from 
other sources, the source of glitches in the HFI bolometers is 
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Fig. 6. Distribution of glitches with respect to the log of their peak amplitude and of ratio of the slow tail to fast part amplitudes. This 
last ratio is estimated by ai/g m - dx , where 5/ is the fitted amplitude of the long glitch template (which has been normalized for long) 
and g max is the glitch maximum. Left and right panels are for one of the PSB-a (143-3a) and one of the PSB-b (143-3b) detectors, 
respectively. The same distributions are shown in linear scale in the upper left comer of the windows. The dash lines indicate the 
observed non-linearity of the amplitude of the slow part of long glitches relative to the peak amplitudes. Long glitches populate 
the distribution along those lines and short glitches populate the distribution along the vertical axis close to zero. Slow glitches 
are clearly visible in the distribution centered at around 6, for the PSB-a, and are absent for PSB-b bolometers. The distribution of 
template amplitudes at the limit of weak glitch peak amplitudes is wider for PSB-b due to larger estimation errors, as long glitches 
for PSB-b (and SWB) have lower amplitude intermediate time constants than for PSB-a. For some bolometer (as for this PSB-a) we 
observe a double peaked distribution of long glitch template amplitudes. 
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Fig. 8. Mean rate of the three categories of glitches long, 
short and slow glitches. These rates are computed using the 52 
bolometers for long and short glitches, the 16 PSB-a bolometers 
for slow glitches. 



from galactic cosmic rays. Indeed, other sources, such as on- 
board radioactivity and solar protons as detected by the mirror 
on WMAP (Jarosik et al. 2007) were suspected to contribute to 
the glitch rate, but would not follow the trends with time that are 
observed. 



2.5. Population counts 

Figure 1 1 shows the distribution dN/dE of the three populations 
of glitches as a function of the amplitude of events in signal to 




Fig. 9. Glitch rate per bolometer. Each point is the mean value 
during a survey (6 months). 



noise units for one of the detectors. The type of each event is 
determined from the output of the method detailed in 4.1. We 
can see that long glitch events are dominant at lower amplitudes 
while short events dominate at higher amplitudes. The separa- 
tion between short and long events is not very efficient below « 
20 sigmas of the noise as the tail of long glitches, used to dis- 
tinguish between the two types of events, is barely detectable. 
Nevertheless, we will see later that long glitches are dominant 
at lower amplitudes. Long glitches are dominant in the overall 
counts, and the tail of those events is responsible of the noise 
excess at frequencies between 0.001 Hz and 0.1 Hz, requiring 
their subtraction in order to reach the expected level of noise 
(see Planck HFI Core Team 2011, and Section 4.2). The distri- 
bution of long glitches is well fitted by a steep power law of 
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Fig. 10. Top : Normalized glitch rate (black line), cosmic rays 
flux measured by the SREM for deposited E > 3 MeV in purple, 
and E > 0.085 MeV in pink, as a function of time. Bottom : 
monthly smoothed number of sun spots. 




Amplitude (signal to noise) 

Fig. 11. Distributions dN/dE of the three families of glitches 
with respect to the peak amplitude in signal to noise units, for 
three different bolometers: 143-la, 143-2a and 217-1. Blue is for 
the short glitch population, black is for long, red is for slow and 
green is for total. We have displayed populations for bolometers 
with very different behaviours: 143-la has one of the least glitch 
rate, 143-2a and 217-1 have a high glitch rate and we observe 
clearly the faint-end break of the glitch counts for 217-1. This 
last bolometer has no slow glitches (see Section 2.3) Power laws 
are shown for comparison in dash lines. Indices and amplitudes 
of power laws are chosen to match the distribution of one of the 
bolometers (143-2a). Indices are -2.4 for long and slow glitches, 
and -1.7 for short glitches. The vertical dash line indicate the de- 
tection threshold. There is no attempt to separate the long and 
short glitch populations below 20 sigmas of the noise. The slow 
glitch spectrum is shown down to 5 sigmas. 



index -2.4 between typical amplitudes of 20 and 1000 sigmas of 
the noise (with some variations of those values with respect to 
bolometers), and show breaks at both bright and faint-ends. The 
faint-end break is very close the glitch detection threshold fixed 
to 3.2 sigmas, and is detected very significantly as we will dis- 
cuss later in Section 3.2. We observe that sub-millimeter chan- 
nels (545 and 857 GHz) are more sensitive to long glitches. For 



those detectors, we observe a peak in the differential counts at 
energies relatively close to the detection threshold. 

We observe more differences between bolometers in the long 
glitch distributions than in the short glitch distributions, as can 
be seen in Figure 18 comparing the cumulative distributions for 
different bolometers and glitch families . The distribution of 
short glitches has a bump at amplitudes of m 3000 sigmas for 
all bolometers from 100 to 353 GHz. There is a clear break at 
very high amplitudes which is mostly due to the non-linearity of 
detectors at those very high amplitudes and to the saturation of 
the analog to digital convertor (ADC, see Planck Collaboration 
VI 2013, for a discussion about the effect on data) for a discus- 
sion about the effect on data). The bump in the distribution is not 
apparent in sub-millimeter detectors, but this is due to the non- 
linearity of detectors reducing it. For amplitudes below w 1000 
sigmas, the distribution of short glitches is well represented by a 
power law with a typical index of -1.7. Short events with an en- 
ergy corresponding to the bump appear to be shorter than events 
corresponding to the power law. The shape of the distributions 
of slow and long glitches are very similar. By computing the 
distributions of the amplitude of the tail of both long and slow 
glitches, the amplitude of both distributions are also very similar. 

2.6. Correlation in time 

The time lapse between two consecutive glitch events exhibits an 
exponential shape as shown in the graph on the left in Figure 12. 
The absolute value of slope fitted between 1 and 2 minutes is 
roughly equal to the computed rate, just obtained by the num- 
ber of events above 10 cr in six months (right part). Clearly, all 
bolometer glitch distributions exhibit a Poisson law. These data 
are from the third survey data when the glitch rate was relatively 
constant and events above 10 cr were selected to avoid any con- 
tamination by false events, solar flares and pile-up. We cannot 
check multiple simultaneous events because we cannot distin- 
guish between a pile-up of multiple events and a single more 
energetic event. Thus the glitch data are compatible with pure 
random events which are not correlated in time and is consistent 
with the cause of events as individual hits by galactic cosmic ray 
particles. 

2. 7. Coincidences between detectors 

We also studied the coincidences of events between bolometers. 
We found that the largest fraction of coincidences, ~ 99%, are 
between PSB-a and PSB-b detector pair in the same module. The 
remaining ~ 1 % of coincidences are particle shower events that 
effect a large fraction of the detectors in the focal plane in dif- 
ferent modules. The PSB pair coincidences are discussed here 
and the shower events are discussed in Section 3.6. At the ar- 
rival time of each glitch detected in a PSB-a bolometer, we have 
measured the level of signal in the PSB-b in the same module. 
The glitches in PSB-a are sorted into the three categories (short, 
long and slow), then we compared the amplitude of the signal in 
PSB-b without distinction on the nature of the signal in PSB-b 
with the amplitude of the glitch at the maximum in PSB-a and 
computed a 2D histogram of the two quantities. The 2-D dis- 
tribution for short, long and slow events in PSB-a are shown in 
Figures 13, 14 and 15, respectively. The axes indicate the log of 
data normalized by the noise rms. The events below 100 sigmas 
are not displayed since the separation between short and long 
glitches is uncertain. 
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Fig. 12. Left: Histogram of time lapse between two consecutive glitch events for each bolometer normalized by total number of 
events. Colors are for different detectors. We observe Poisson laws as expected for random events . Right: glitch rate derived from 
slope of time lapse distribution compared to the mission averaged actual rate for each bolometer. 




Fig. 13. Representative short glitch coincidence in a PSB pair. 
The x-axis represents the log of the amplitude relative to the 
noise of detected glitch in the PSB-a bolometer, the y-axis is the 
log of the value of the signal to noise in the PSB-b bolometer. 
The grey scale indicates the log of the number of events in log- 
arithmic bins. Detected events are not shown below 100 sigmas 
in PSB-a, as we expect non-negligible contribution from long 
glitches, the separation between long and short being uncertain 
at low amplitudes. A large fraction of events are in coincidence 
between PSB-a and PSB-b. 



Fig. 14. Representative long glitch coincidence in a PSB pair. 
See description in Figure 13 caption. Events above 50 sigmas are 
selected such that they match the long glitch template. All events 
without selections are taken below 50 sigmas, as the long glitch 
population dominate at low amplitudes. The solid line indicates 
identical signal in both PSB-a and PSB-b. Most (if not all) long 
events as associated with significant events in PSB-b, and a large 
fraction have approximately the same amplitude in PSB-a and 
PSB-b. 



We observe that most of the long glitches detected in PSB-a 
are also detected in coincidences in PSB-b. Furthermore, we find 
that for amplitudes above 300 sigmas, all glitches are detected 
in coincidence, as all events are seen in PSB-b above « 0.3% of 
the amplitude in PSB-a (or reversely). We cannot state this for 
events below 300 sigmas since a significant fraction of the coun- 
terparts in PSB-b are below the noise level. Nevertheless, the re- 



sults of the coincidence analysis are compatible with the fact that 
all long events are in coincidence between PSB-a and PSB-b. A 
large fraction of these glitches is seen with an amplitude which 
is proportional in the two detectors, as those events lie along a 
line in the diagram with some intrinsic scatter, meaning that a 
large fraction of long glitches are coincident with the same ratio 
of amplitudes between 'a' and 'b'. The factor of proportionality 
is almost 1 in signal to noise units for this PSB pair, but can vary 
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Fig. 15. Representative slow glitch coincidence in a PSB pair. 
See description in Figure 13 caption. No event is shown below 
50 sigmas to avoid leakage from the long population. All slow 
glitches in PSB-a above * 10 3 sigmas are associated with sig- 
nificant events in PSB-b. Below « 10 3 sigmas the majority of 
events are associated with significant events in PSB-b, for the 
remaining fraction the fact that we do not see a significant signal 
might be a threshold effect. 



between 0.25 and 4 depending on the pair (see Section 3.2 for 
a discussion). A very significant fraction of long glitches have 
counterpart in the other associated bolometer, but with a differ- 
ent amplitude, which is lower most of the time. The number of 
events in the distribution along the line of one to one correlation 
is of the order of 50% of the total number for events with am- 
plitudes lower than * 300 sigmas of for events above * 1500 
sigmas, and is of the order 35% for ampltitudes in between. 

The distribution of PSB-a versus PSB-b for short glitches is 
very different than the expected distribution for random coinci- 
dence. The expected distribution of random coincidences is eval- 
uated by measuring the statistics of the signal in PSB-b 10000 
samples (about one minute) after each event in PSB-a, and is 
shown in Figure 16. In particular, we can see a wide cloud of 
points in the 2D histogram above 10 sigmas in PSB-b. A frac- 
tion of 50 % of the events populate the region of random coinci- 
dence. This is also shown by Figure 17 displaying the distribu- 
tion of the signal in PSB-b for events in PSB-a with amplitudes 
between 300 and 1000 sigmas. So there are strong indications 
that about 50 % of short glitches are seen in coincidence between 
PSB-a and b. For a given amplitude of short glitch in PSB-a, we 
observe a wide distribution of amplitudes in PSB-b. 

All very high short amplitude glitches have significant coun- 
terparts in the other bolometer at the level of * 1/1000. This 
is attributed to the cross-talk between the two bolometers. This 
will be discussed in Section 3.4. 

As already stated, slow glitches are seen only in PSB-a. 
Nevertheless, there is a significant signal in PSB-b associated 
with slow glitches in PSB-a. Indeed, we can see that the distri- 
bution of counterpart in PSB-b is very different than the random 
distribution shown in Figure 16, which means that slow events 



Fig. 16. Expected distribution for random coincidences. See de- 
scription in Figure 13 caption. The value of the signal to noise 
in PSB-b for each short event in PSB-a is taken 10000 samples 
after. We use short events in PSB-a such that this distribution can 
be directly compared to the coincidence of short glitches. 




Signal to noise 143-4b 



Fig. 17. Distributions of the signal values in PSB-b for short 
events with peak amplitudes between 300 and 1000 sigmas in 
PSB-a in red. Random distribution of values in PBS-b in black. 
The excess of the red curve over the black for amplitudes above 
3 sigmas represents about 50% of the total area. Thus, at least 
50% percent of the events are in coincidence between PSB-a 
and PSB-b. 



in PSB-a are associated with signal in PSB-b. The counterparts 
generally have a wide range of amplitudes for a given energy of 
slow glitches which makes difficult the measurement of associ- 
ated time constants. Visual inspection of those events allow to 
say that there is a variety of counterparts, some have faster de- 
cay than slow glitches and have similar rising time of m 10 ms, 
some others (the larger energy ones) are similar to long glitches. 
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Fig. 18. Cumulative counts (in hour" 1 ) for the three categories 
of glitches and for several bolometers. Blue curves are for short 
glitches after non-linearity correction above typically 10 keV (as 
described in the text); black curves are for long glitches; red 
curves are for slow glitches; and green curves are for the total 
counts without non-linearity correction, which is why they lie 
below the blue curves. Glitches are not separated into the dif- 
ferent categories below 10 sigmas of the noise (« 20 eV). All 
curves have been recalibrated in energy for each bolometer such 
that the distributions of short glitches match at energies around 
lOkeV (we use the same coefficient for all populations), the ab- 
solute calibration being fixed for one arbitrary detector. Notice 
that the distributions of short glitch match very well. There is a 
large scatter in the energy scaling for the other populations of 
glitches. 



3. Interpretation 

For a complete understanding of the long and short glitches ob- 
served in the HFI flight data, several laboratory tests have been 
performed using 23 MeV protons from a TANDEM accelerator 
and two radioactive alpha sources. The detailed results of all the 
performed tests are presented in the HFI special paper Catalano 
etal. (2013). 

3.1. Interaction with the grid/thermistors 

There is clear evidence that the short events are resulting from 
cosmic rays hitting the grid or the thermistor. Indeed, those 
events have a fast rising time and have a fast decay, and the trans- 
fer function built from the short glitch template (see Catalano et 
al. 2013, for a comparison), is in good agreement with the HFI 
optical transfer function (Planck Collaboration VII 2013), so the 
energy must be deposited in the environment close to the ther- 
mistor. Figure 18 shows the cumulative counts, i.e. the number of 
events N(>E) with an amplitude higher than a given value E per 
hour, of short glitches after converting the maximum amplitude 
of each glitch in units of energy (keV) for one reference bolome- 
ter (100-lb). The energy calibration for the other bolometers has 
been determined such that the cumulative distributions match the 
one for the reference bolometer at energies of about lOkeV. For 
the conversion we use the measured heat capacity of the NTD 
crystal + grid system, as we assume that the energy of the events 
is deposited in this system. The relative fitting of the energy is 
necessary for the comparison between bolometers because of un- 
certainties in the heat capacity of the system and also in how 
much the peak amplitude is reduced after averaging in the inter- 



val of one sample with the electronics. Faster bolometers have 
glitches with more reduced amplitudes than slower bolometer. 
The effect is important for sub-millimeter channels for which 
the time constant is very small. We evaluated that the relative 
correction factor should be of the order of three between 100 
and 545 GHz bolometers. We have also reconstructed and cor- 
rected the counts for high amplitudes which are affected by non- 
linearity (typically above w lOkeV) by readjusting the values of 
event energies from the measurements of the slow tails of short 
events, which are not affected by non-linearity, instead of tak- 
ing directly the peak amplitude. This allows to reconstruct the 
distribution more than a factor of ten above the saturation by 
the electronics. We can see that the counts of short events for 
all detectors match with a very good accuracy, around the en- 
ergy bump in particular, considering the fact that no rescaling is 
performed on the y-axis, but only on the energy axis. 

Events in the bump are expected to be associated with cosmic 
rays hitting the thermistor. Since thermistors are identical for all 
bolometers, we expect the same rates for all bolometers, this is 
precisely what we observe, as shown in Figure 18. Given the di- 
mensions of the thermistors of 30 x 100 x 300 /im 3 , and the fact 
that typical galactic protons of energies around 1 Gev (peak of 
the spectrum) deposit the minimum energy of a: 1 .47 MeV cm 2 /g 
in the Germanium, which is the minimum of the stopping power 
for normal incidence in the largest surface, we expect typical de- 
position energy of w 15 keV in the thermistor (the first energy be- 
ing dominant since it corresponds to the largest cross-section). 
This is also in agreement with the energy measured experimen- 
tally, given the fact that the energy calibration is uncertain within 
~ 50%. The amplitude of the bump in the cumulative distribu- 
tion which is of about 10 events per hour, also matches with 
the expectations. Indeed, by integrating the expected cosmic ray 
spectrum at L2, and considering that low energy particles are ab- 
sorbed by the spacecraft, we can predict a total number of par- 
ticles of N ~ 5 particles/min/mm 2 . This leads to 0.1 glitch/min 
on the thermistor by taking the largest surface. 

The events depositing an energy of about 10 eV to 1 keV and 
populating the power law as seen in Figure 18, are expected to 
be the result of cosmic ray hitting the grid. Calculations predict 
a total number of events on the grid from 0.4 to 1 .2 glitch/min, 
given that the surface from the grid varies from 4.77 • 10 5 fxm 2 to 
3.95 ■ 10 4 fim 2 depending on the bolometer. We observe a higher 
number of events by a factor of 2 to 4, but as well we see in 
the next paragraph, about half of the events for PSB bolometers 
might result from electrons ejected from the other grid by the 
impact of a cosmic ray proton. We observe some dispersion in 
the amplitudes of the power law distributions from bolometer to 
bolometer which is not obviously correlated with the area of the 
grids. This is attributed to the uncertainties in the intercalibration 
of energy to match the counts, as detailed at the beginning of this 
section. 

The coincidences analysis, described in Section 2.7, of short 
glitches (Figure 13) at intermediate amplitudes which populate 
the power law, show that a fraction very close to 50% of events 
are seen in both PSB-a and PSB-b with no phase shift. A possi- 
ble model is that events hitting one of the grid are ejecting some 
electrons hitting the other grid. The rate of 50% correlation can 
be explained by geometrical effects like particles coming from 
the top or bottom; a particle coming from the top could hit the 
PSB-a grid, extract electrons and project those onto the PSB- 
b grid, whereas a particle coming from the bottom and hitting 
the PSB-a grid would not project electrons on the PSB-b. The 
observations indicate then that nearly all events hitting one of 
the grid and coming from the appropriate direction must eject 
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electrons hitting the second grid and depositing about 10% of 
the energy deposited in the first grid. We also observe that some 
events are hitting both grids. Those appears in the cloud of points 
with the same amplitude in both PSB-a and b. Those represent 
about 2% of the events. We observe a lower coincidence rate for 
events at higher amplitudes corresponding to the bump. This is 
explained with the hypothesis that those events correspond to di- 
rect impact on the thermistor, which is not aligned with the grid 
of the other bolometer. Nevertheless, we have seen that all those 
high energy events are in coincidences with a small amplitude 
event in the second bolometer. This is attributed to the cross- 
talk between bolometer in a pair, as discussed in Section 3.4. 
The non-linearity, appearing with the curved shape in the coin- 
cidences, can be entirely explained by the saturation of the high- 
est amplitude glitches. Nevertheless, we can evaluate the level of 
cross-talk as between 0.01 to 0.44% depending of the bolometer 
pair, for « 3000 sigmas events for which the effect of saturation 
is small. The high amplitude events correspond to an energy of 
as 15keV left on the grid+thermistor by a particle, this is very 
close to the expected energy deposited by a 1 GeV proton on the 
thermistor. 

3.2. Interaction with the wafer 

We identify the long glitches as produced by cosmic rays hit- 
ting the silicon die. This was first indicated by the ground tests 
(Catalano et al. 2013), showing that the NTD thermometer is 
sensitive to a temperature change of the silicon die. The HFI 
ground-based calibration show a rate of events compatible with 
the cosmic rays flux at sea level over the silicon die surface and 
also that almost all these events are in coincidence between PSB- 
a and PSB-b. The understanding is the following: phonons gen- 
erated by the event impact in the silicon die and produce fast ris- 
ing time of the Germanium temperature, which decays with the 
bolometer time constant. The slow part is the thermal response 
of the entire silicon die temperature rising and then falling as 
the heat conducts out from the die to the heat sink (Catalano 
et al. 2013). This hypothesis is reinforced by the comparison 
of cumulative counts N(>E) of long glitches from bolometer 
to bolometer which is shown in Figure 19. To make this com- 
parison, we have normalized the counts such that we evalu- 
ated the number of event per unit of time and per unit surface 
of the silicon wafer (as the total surface varies from bolome- 
ter to bolometer). Also, as for short glitches, we have intercal- 
ibrated the energies of events between bolometers by match- 
ing the different counts at measured energies around 0.05 keV, 
which correspond to the deposited energy on the silicon die of 
« 10 3 keV after absolute calibration. This absolute calibration 
is performed such that the observed faint-end break of counts, 
which is clearly observed in the figure, match the expected min- 
imum deposited energy by ~GeV protons on the silicon die, 
which is * 140keV, for normal incidence, as described later. 
Energy calibration cannot be performed accurately from mea- 
surement due to important uncertainties of heat capacities of the 
silicon die. Figure 19 shows the cumulative counts for all HFI 
bolometers used for the scientific analysis. We clearly see that 
they all match in shape and amplitude with very small scatter. 
For energies below 20 sigmas of the noise, we have also com- 
puted the cumulative counts for all events without distinction on 
families (green curves in the figure), but we have already seen 
that long events dominate at low energy, so the total counts are 
representative of counts of long glitches. The faint-end break in 
the counts is detected without ambiguity. This is an important 
result since it shows that there is a limited number of low en- 




Amplitude (KeV) 

Fig. 19. Cumulative distribution N(>E) of long glitches per 
unit surface of the wafer for all bolometers used for the scien- 
tific analysis. The black solid lines correspond to selected long 
glitches. The green lines correspond to total glitches, but are rep- 
resentative of long glitches since they dominate at low energy. 
A relative calibration of the energy, which indicated in the x- 
axis, has been for all bolometers such that the cumulative counts 
match at a; 10 3 keV. The absolute calibration of the energy is 
performed such that the faint-end break of the counts match the 
expected deposited energy of * 140 keV in the silicon die. This 
calibration on the model is necessary since we only measure 
the energy in the grid+thermistor system, and the heat capacity 
and link conductivities are not known with a sufficient accuracy. 
Also, we use the values of peak of glitches to record the glitch 
amplitudes which are not representative of thermal processes, as 
the fast part of long glitches are expected to result from ballis- 
tic phonon effects. Nevertheless, we observe that the shape of 
the cumulative counts match for all detectors with a low disper- 
sion. The break of the counts at faint-end is very significantly 
detected. This is our best indication from flight data that long 
glitches result from particles hitting the bolometer silicon die. 
Dot dash and dash lines correspond to the predicted power law 
distributions of deposited energy from protons, and helium, re- 
spectively, from the toy model developed in 3.3 for sensibly en- 
ergies higher than the stopping power energies. Predicted min- 
imum deposited energies for helium and protons are indicated 
with doted lines. We interpret the second apparent bump at en- 
ergies ~ 20 times the energy of the break as the signature of 
He nuclei. The fact that we observe an excess in the counts at 
the expected energy of the minimum stopping power for helium 
reinforces our hypothesis for the origin of long glitches. 



ergy glitches in data which are undetected, and would affect the 
cosmological results. From this limit, we measure that the total 
number of events penetrating the shielding of the satellite around 
the bolometer is about 4.5 per second and square centimetre. 
This is in very good agreement with the expected value in the 
bolometer environment computed in Section 3.1 of « 5 s^'cm" 2 . 
The slope and amplitude of the distribution can be predicted with 
a simple model of interaction of particle with the silicon die as 
shown by the toy model described in Section 3.3. Predictions are 
also represented in Figure 19 for primary galactic protons and 
Helium nuclei. We observe a very good match of the model with 
the data showing that particles detected by Planck at the ener- 
gies above the faint-end break are primary galactic protons. The 
barely apparent second bump in cumulative counts at deposited 
energies ~ 3000 keV probably correspond to the contribution 
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from Helium nuclei, which are expected to contribute to ss 10% 
of the counts at those energies. High energy galactic electrons 
could also contribute to the measured counts for high deposited 
energies 

The very large scatter in the distributions of long glitches 
with respect to different detectors, calibrating the energy on the 
short distributions, as shown in Figure 18, is due to the variation 
of thermal links between the grid and the silicon die, as well as 
silicon die heat capacity, reducing or increasing the fraction of 
energy transmitted to the thermistor for both phonons and ther- 
mal processes. The detectors with smaller glitch counts (e.g 143- 
la and 143-5 in Figure 8) are those with smaller fraction of the 
energy propagating from the wafer to the grid. For those detec- 
tors we do not see the break in the counts which should be below 
the detection threshold. We observed that in general PSB-b have 
are more sensitive to long glitches than PSB-a (a few exceptions 
exist). 

Detectors with smaller rates of glitches should have more 
undetected glitches than detectors with a high rate. This hypoth- 
esis is verified with the help of the cross-correlation analysis be- 
tween bolometer signal from the same PSB pair described in 4.2 
(see Figure 27 in particular). We observed that bolometers pairs 
with the smaller rate of glitches are the ones with higher corre- 
lations around 1 Hz, which is attributed to glitches below thresh- 
old. Whereas pairs for which one of the bolometer has a visible 
faint-end break in the counts has smaller correlation of noise. 
The impact on final data after processing is studied in 5.2. 

We have already seen that events hitting the bolometer sili- 
con dies must deposit a faction of their energy, which is typically 
the minimum stopping power of the silicon for events close to 
the normal incidence and with energies of 140 keV. Then most 
of the events must hit instantaneously the two wafers in PSB 
pairs. This is compatible with what we observe in the coinci- 
dence graph shown in Figure 14, Section 2.7. Indeed, we find 
that between a third and half of the events are depositing the 
same energy in both PSB-a and PSB-b with some intrinsic dis- 
persion. We mentioned in 2.7 that the factor of amplitude be- 
tween events in signal to noise units varies from 0.25 to 4 de- 
pending on the bolometer pair, but this factor becomes 1 after 
recalibrating in units of deposition energy on the silicon die as 
expected. Furthermore, 100% of the events (at least for ampli- 
tudes higher that 300 sigmas) happen in coincidence between 
the two bolometer pairs but with different energy deposition, by 
a factor of 300 at maximum. This difference could be explained 
by the fact that some relatively low energy events loose a sig- 
nificant fraction of their energy in the first silicon die and then 
deposit more in the second, as the stopping power is a decreasing 
function of the energy for the considered particles. In principle 
we should have a small proportion of events hitting the corner 
of the bolometer wafer with no counterparts in other bolome- 
ters, but we do not have evidences of such events. A complete 
modeling of the interaction of particles with the two silicon dies, 
including accurate physical modeling of the interaction is post- 
poned to a forthcoming publication. 

Very importantly, the coincidence study indicates that the 
contribution to the glitches by secondary particles, like e.g., delta 
electrons is negligible since those low energy events would de- 
posit all their energy in the first wafer without counterparts in 
the second. 

The slow part of long glitches is related to the thermal relax- 
ation of the energy. The different time constant can in principle 
be related to the different components of the bolometer. This will 
be studied in the future publication Spencer (2013). 




Fig. 20. Geometrical scheme of the toy model of cosmic ray in- 
teraction with the silicon die. 

3.3. Toy model for interaction with the silicon die 

We have developed a toy model considering the impact of cos- 
mic rays at the second Lagrange point with the silicon die. We 
start with a solid square box made of silicon with the same equiv- 
alent surface of the real silicon die (see Figure 20). We consider 
that the side of the square is much greater then the thickness of 
the silicon die. The input of the model are: 

- Geometrical parameters of the bolometers. 

- Stopping power function and density of the silicon die. 

- Energy distribution of CRs at L2 (see e.g., Adriana 201 1). 

By integrating over the solid angle, the surface and the inte- 
gration time, we obtain an analytical equation as: 
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where No is the amplitude of the spectrum of incoming protons 
or alpha particles at L2, At the integration time, E p +„ the refer- 
ence proton energy, y and j3 are the power law indexes of the fit 
between the stopping power function and the energy distribution 
of the proton at L2 respectively, Eo = p s u ■ d ■ SP(E p + ) is the 
reference absorbed energy for an orthogonal impact, p„y the den- 
sity of the silicon, d the thickness of the silicon die, S P(E p + ) the 
stopping power function calculated at a reference proton energy, 
and E^ie is the energy absorbed by the silicon die. 

The predictions of this model for protons and alpha particles 
are shown in Figure 19, together with the energy spectrum mea- 
sured in bolometer data. The power law index fit well with the 
flight data and the model is able to cover almost all the range 
of energies. We conclude that in terms of rate and energy dis- 
tribution the thermal coupling between the silicon die and the 
NTD thermometer likely explain the long glitches seen in HFI 
in-flight data. 

3.4. Cross-talk between bolometers 

Very high amplitude short events in one of the bolometers from 
a PSB pair is always seen in coincidence with a small event 
in the other bolometer. This is attributed to cross-talk between 
bolometers. The estimated level is « 10~ 3 as can be seen in 
Figure 13 for the considered bolometer pair, and varies be- 
tween 0.01 to 0.44% depending of the pair. Those events show 
a phase shift which varies from 0.25 to 1.4 samples depending 
on the pair, and the transfer function of the cross-talk signal is 
different than the primary glitch transfer function. We observe 
longer time-constants than the bolometer for the cross-talk sig- 
nal, which scale from 10 to 30 ms. The origin of cross-talk be- 
tween detectors from the same PSB pair might be electrical or 
thermal. There is no clear evidence at this point. 
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Fig. 21. output of SREM diodes TCI, TC2 and TC3 (left axis) 
and temperature control heater power on the 0.1K plate (right 
axis) and as a function of time for large solar flares on March 7 
(top) and August 4 (bottom) in 201 1 . 



We have attempted to estimate cross-talk signal between 
bolometers which are not mounted in the same pair. To do so 
we have averaged the bolometer signal at locations of detected 
high energy events in one of the bolometers. We found a signif- 
icant but very small cross-talk signal of the order of 3 x 10~ 5 
for some couples of detectors within the same electronics belt. 
However, we did not measure any significant signals in data from 
other bolometers, and we measured a higher limit of 3 x 10 6 . 

3.5. Solar Flares 

For most of the mission, solar activity was remarkably low 
(Mewaldt et al. 2010). In 2011, for the first time in the mis- 
sion, there were several large solar flares. These flares provided 
a useful test signal to correlate the signal measured on the out- 
side of the spacecraft with the SREM to signals due to particle 
impact on HFI. The glitch rate noticeably increases during so- 
lar flares. The heater power used to regulate the 0.1K plate de- 
creases with increasing signal, and hence particle flux, measured 
by the SREM diode sensors as shown in as shown in Figure 21. 
The time scale of each flare was slow enough that the temper- 
ature control loop was able to compensate for the bulk heating 
produced by the increased particle flux as evidenced by the small 
phase shift of the heater power compared to the onset of signal in 
the SREM. Note in Figure 21 that the peak measured signal for 
TCI and TC3 are only ~ 40% for the two flares. The signals of 
the heater power and TC2 are similar to each other in each flare. 
However, the peak signal of each is very different comparing the 
two flares. In addition, there is structure in the signal for TCI and 
TC3 that are not in TC2 or the heater power response. We find 
that this correlation holds between the heater power and TC2 for 
all flares. The diode TC2 has the most shielding of the 3 diodes 
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Fig. 22. Coincidence histogram (solid line) processed for the en- 
tire mission, based on the result of the glitch detection algorithm. 
Random coincidences dominate for low number of bolometer 
touched. Coincidences above 15 bolometer touched are real cor- 
relations. 



in the SREM, 1.7 mm of Aluminum and 0.7 mm of Tantalum 
which passes only ions and protons with energies > 39 MeV. The 
other diodes are shielded by 0.7 mm of Aluminum for TC3 and 
1 .7 mm of Aluminum for TC 1 . This demonstrates that the space- 
craft and instrument surrounding the bolometers shield particles 
with energies at least up to 39 MeV and all solar electrons. This 
is equivalent to the stopping power of ~ 1 .5 cm of Aluminum 
(Bergeret al. 2012). 

3.6. High Coincidence Events 

Since the beginning of the Planck survey, high temperature 
increases of the 0. 1 K plate have been observed through the 
bolometers and thermometers. These events have two specific 
characteristics. There is a high number of glitches detected in 
the same time in multiple bolometers and there is an increase of 
the focal plane temperature of nK to fiK. In this section we will 
analyze the characteristics of these High Coincidence Events 
(HCE). 

3.6.1. Detection 

HCE are identified by detecting a large number of precur- 
sor glitches in coincidence in different detectors (so called 
"touched" bolometers) on the 0.1K plate. These coincidences 
are a precursor to a temperature excursion of the 0.1 K plate. 
Other automated techniques that focused on finding only the 
temperature excursion were applied to the data, but were not 
successful at identifying HCE in the time stream. One reason 
for this is that for most of the events, this temperature rise is 
below the thermometry noise level. 

Figure 22 shows the histogram of the number of events in 
coincidence for the full mission for 15 ms bins. The distribution 
at low number of coincidences is nearly compatible with random 
coincidence given the measured rate of glitches per bolometer. 
The distribution deviates from random coincidences above 13 
bolometer touched. We set the detection threshold at 15 in order 
to select only real coincidences linked with HCE. 

We detect * 100000 HCE for a number of coincidences 
greater than 15 or on average 5 events per hour. 
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3.6.2. The Different Types of HCE 

For each event above the thermometry noise level we calcu- 
late separately the median of the temperature time stream for 
"touched" bolometers and bolometers without a precursor glitch 
or "untouched" bolometers. We can separate events in two cate- 
gories : fast and slow, which differ in the rise time of the temper- 
ature shown in Figure 23. Depending on the type of the events, 
the maximum temperature value is calculated between 8-12 and 
80-120 seconds after the coincidence respectively for fast and 
slow events. 

For fast HCE, coincident glitches occur in several, but not 
all bolometers and are grouped together on the 0.1K plate. The 
temperature rise time constant is a: 5 seconds. The decay time 
constant is w 20 minutes and matches with the thermal time con- 
stant of the 0.1 K plate. 

The precursor glitches for a fast HCE are short or long 
glitches. The separated "touched" and "untouched" bolometer 
temperature time stream stacking average for fast HCE greater 
than 1.5 yuK is shown in Figure 23. The "untouched" curve rises 
slowly, but decays in phase with the "touched" after all bolome- 
ters thermalize with the entire 0.1 K plate. The distribution of 
fitted temperature amplitude for all the fast HCE (coincidences 
greater than 15 bolometers) obtained by template fitting is given 
in Figure 24 left panel. We observe that fast events above 0.8 /iK 
are individualy detected. For events detected with signal to noise 
less than one, we find a non-zero amplitude of 0.057 muK. We 
also observed that the rate of fast HCE is modulated as cosmic 
rays over the mission. 

The facts that the "touched" bolometers are grouped in the 
focal plane and precursor glitches of fast HCE are the same as 
for individual particle hits directly on the Si die (long) or grid 
(short), we conclude that these events are due to showers of sec- 
ondary particles over part of the focal plane. The temperature in- 
crease of the 0.1K plate results from the low energy secondary 
particles from the shower stopped by the bolometer plate. The 
rate of the most energetic events is consistent with the geome- 
try of the 0. 1 K plate and flux of Galactic cosmic rays energetic 
enough, » GeV, to produce such showers (Adriana 201 1). 

For slow HCE, the coincident precursor glitches occur in 
all, or nearly all, bolometers. As shown in Figure 23, there is 
a small temperature decrease of » 1 second after the precursor 
glitch and before the temperature rise with a time constant of « 
30 seconds and, as observed for fast HCE, a decay time constant 
of * 20 minutes. This initial decrease prior to the rise is puzzling. 

The precursor glitches for these HCE are very specific and 
seems to be different to any categories studied in the Section 2. 
They have a rise time of ss 10 ms and a decay time which varies 
from event to event but of » 100 ms which is a lower value from 
what is seen for the slow glitches (see Section 2.3). The dis- 
tribution of fitted temperature amplitude for all the slow HCE 
is shown in Figure 24 (right panel) and the validation thresh- 
old value is 1.3//K. Slow HCE have a less steep distribution 
than fast HCE and dominate above 2//K. As opposed to fast 
HCE, we have observed that the number of slow events de- 
crease along the mission. We use the measured thermal prop- 
erties of 0.1 K plate(Planck Collaboration II 2011) to estimate 
that slow HCE correspond to energy deposition of « 1 TeV with 
a rate which is not consistent with the galactic cosmic ray spec- 
trum( Adriana 2011). We have considered elastic relaxation of 
cracks and anomalous response of the 0.1K stage PID loop as 
possible explanations. However, we have yet to reproduce a slow 
HCE experimentally on the flight focal plane, for example with 
temperature or heater steps, or in analogous ground tests. At this 



time, we have not assigned a physical cause for the slow HCE 
events. 



3.7. Particle interaction with bolometer plate 

It has been shown in Planck Collaboration II (201 1) that the in- 
teraction of particles with the bolometer plate are responsible 
of the correlated low frequency noise between bolometers. This 
excess noise was not seen in ground measurements (Pajot et al. 
2010) and is modulated with solar activity the same as the glitch 
rates. 

Although the mean temperatures of the 0.1K and dilution 
cooler plates are stabilized with great precision, it is only at a 
single point near the control heater. Distributed heating from par- 
ticles and other sources cause gradients within the plate. As de- 
scribed in Planck HFI Core Team (201 1), the signal from the 2 
dark bolometers was used to track the temperature of the 0.1K 
plate, independent of the control thermometer. The signal from 
each was smoothed over a running window of 2 minutes (or ~ 2 
full rotations of the spacecraft) and combined to form a thermal 
template which is then cross correlated with each bolometer and 
used to remove long term drifts from the bolometer time stream. 
Both the control heater power and the thermal templates are cor- 
related with diode sensors in the SREM indicating a significant 
source of heating (although small <10nW) was due to heating 
from interaction with the particle flux. 

We have evaluated the contribution of fast HCE to the 
bolometer noise power spectrum. We start by assuming that the 
distribution of the amplitudes of fast events follows the fitted 
power law shown in Figure 24. Since the power law has a high 
index (« 3), some break have to be put at the low amplitude 
limit. By choosing a limit about 10 times lower than the am- 
plitude such that the total counts match the detected number of 
coincidences for 15 bolometer touched, we found that fast HCE 
might account for approximately 10% of the correlated low fre- 
quency noise between bolometers. By selecting events above 15 
coincidences only, we are not counting all actual coincidences 
since coincidence of a few bolometers, < 10, are overwhelmed 
by the random events. We have tried techniques, such as fitting 
HCE templates over a sliding window of 0.1K plate tempera- 
ture measured by the thermal template, but have not found a one 
successful at detecting the bulk thermal response of the 0.1K 
plate from fast HCE without the precursor coincidence flagging. 
However, just a modest increase in the total HCE events with 
< 10 bolometers in coincidence could explain all of the corre- 
lated noise. We are pursuring models and analysis techniques to 
quantify more accurately the number and effect of these unde- 
tected coincidences. 



3.8. Other events 

The origin of slow glitches is not entirely understood at this time. 
Lab tests using alpha particles were unsuccessful at reproducing 
slow glitches. Nevertheless, the similarities with long glitches 
suggest that both populations have the same physical origin, and 
that they correspond to some energy injected in the wafer. A key 
feature of the slow glitch is that they only occur in the PSB-a. 
The mechanical mounting of the PSB-a bolometer is different 
than the PSB-b bolometer and the SWB bolometer. Lab tests 
using heaters and thermometers mounted at different locations 
near the PSB-a are in progress to attempt to identify any unique 
thermal feature that could cause a slow glitch. 
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Fig. 23. Stacking of 69 fast HCE for bolometers detected as "touched" and "untouched" and 72 slow HCE for all bolometers. 
For fast HCE (left panel), the limit is clear, and the stacking of "untouched" bolometers gives the thermal time constant of the 
copper/stainless steel bolometers plate heating. For slow HCE (right panel) we observe a temperature decrement lower than 1 fjK 
during * 1 second before the temperature rising which is still not understood 
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Fig. 24. Histogram of amplitudes of fast (left panel) and slow (right panel) HCE. Amplitudes are computed by joint fitting of fast 
and slow templates. Solid lines correspond to power laws fitted on the tails of distributions (above 1 yuK), the high amplitude fast 
events above 7 fiK being removed from the fit. The number of those high amplitude fast events is not compatible with the power law 
distribution. Dot dashed lines are gaussians fitted on the distributions near zero amplitude and are representative of fitting errors. 
Parameters of the fitted distributions are given in the figure legend. Fast events (and to some extents slow HCE) have a varying 
temporal shape from event to event. It appears that the variation can be fitted by a contribution of the slow HCE event template. This 
explains in part the negative tail of the distributions of fitted amplitudes. The other reason is the confusion between events inducing 
an anti-correlation of amplitude parameters. 



One of the bolometers (100- la) shows another category of 
glitches which is not detected in any other bolometers. Those 
particular events have a slow response which is similar to the 
slow tail of short glitches, but do not show a fast decaying part. 
They also have a rising time-constant of w 10ms. Those events 
are not accounted in our glitch subtraction method and this in- 
duces a small excess in the residual power spectra at frequencies 
around 0. 1 Hz (see Figure 26, and Section 4.2 for a discussion of 
the impact of such events). The excess is accounted in the total 
noise budget for cosmological studies (Planck Collaboration XV 
2013). 

Finally, precursors of slow HCE are events which are 
only seen in coincidence and associated with slow HCE (see 
Section 3.6 for a description). They are not accounted by the 
glitch subtraction method but are rare enough such that they do 
not affect data analysis. 



4. Glitch Detection and Removal 

A previous version of the method of detection and removal of 
glitches was described in detail in Planck HFI Core Team (201 1). 
Here, we summarize the method and note the changes that have 
been implemented. We also show evaluation of the efficiency of 
the method and the impact on the noise power spectrum in the 
timestream. 

4.1. Method 

The method iterates between sky signal estimation at the ring 
level and glitch detection and fitting. Data are treated indepen- 
dently ring by ring (see Planck Collaboration 1 201 1, for the def- 
inition of a ring). The simple digital three-point filter (0.25, 05, 
0.25) is used to reduce ringing around large events resulting from 
the electronic filter (Lamarre et al. 2010) and demodulation of 
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data. The three-point filter also enhances the efficiency of glitch 
detection since the bolometer time constant is very close to the 
sampling period. 

4.1.1. Glitch detection and template subtraction 

After estimating the sky signal as described in part 4.1.3 and re- 
moving it from the data, events above a noise threshold set at 
3.2 sigmas on a sliding window of a thousand samples are de- 
tected. 2 On each window, a joint fit of the amplitudes of short 
and long glitch templates is performed simultaneously for all 
detected events in the window, centered at the maximum event. 
The short and long glitch templates are built independently of the 
method, and are each a sum of exponentials, three terms for short 
and four terms for long glitches. Templates were estimated using 
events detected by the method used at previous iterations of the 
data processing. We observed that parameters of the templates 
are essentially independent of the changes on the method since 
templates are determined using the list of detected high energy 
glitches which is not sensitive to the details of the processing. 
Consequently, iterations on the glitch template estimation was 
unnecessary. The method used to estimate the templates is de- 
scribed after presenting the main steps of the method, in 4. 1 .4. 
The joint template fitting of both short and long glitches is a main 
improvement as compared to the method published previously 
(Planck HFI Core Team 2011). It has improved the accuracy of 
determination of the long glitch template amplitude and led to 
smaller residuals from glitch in the final noise. Only one overall 
amplitude is used as free parameter for each glitch type and for 
each event. 

Detected glitches are fit to the templates starting from three 
samples after the maximum, for small amplitude glitches in 
PSB-b or SWB, and is increased to eight samples after the maxi- 
mum, for large amplitude glitches, and from six to eight samples 
in PSB-a (because of the presence of slow glitches). In order 
to break degeneracies between template amplitudes in case of 
high confusion between events, we use the prior that amplitudes 
are positive. The parameter fitting is performed with the ana- 
lytical minimization of the x 1 as detailed in Planck HFI Core 
Team (2011), Equation 2. Priors are used by fixing parameters 
for which the values obtained by minimizing the x 2 is outside 
the authorized range. Those parameters are fixed to in particu- 
lar if a negative amplitude is obtained and another minimization 
is performed. 

Slow glitches are detected as long glitches with an extremely 
large overall amplitude. Consequently, those events are pro- 
cessed by the method in the same way than long glitches with a 
higher fitted amplitude relative to the peak. When a slow event is 
detected, another iteration of the x 1 minimization is performed, 
for which the long template is adapted for this event by removing 
the fastest exponential. 

If the measured long glitch template amplitude for each event 
in the vicinity of the largest event is above the threshold, fixed 
to 0.5 times the expected amplitude for a long glitch at the mea- 
sured peak amplitude, then the long glitch template is removed 
from the data. The expected amplitude of the long tail of glitches 
account for the non-linearity described in Section 2.2. To do so 
we use a simple empirical quadratic law which is fitted to the 
data, relating the glitch template amplitude and the peak am- 
plitude. To remove the template from the data, we use the fit- 



2 The choice of 3.2 sigmas is a compromise between false event 
detection and glitch detection completeness as described previously 
(Planck HFI Core Team 201 1), see Figure 7 in that reference. 



ted value of the amplitude of the long glitch template for events 
above 10 sigmas. Glitches with overall amplitude below 10 sig- 
mas are treated either as long or slow glitches, and the amplitude 
is fixed to the expected values since 1.: the long glitches dom- 
inate at low amplitude, and 2.: the fitting errors are larger than 
the model uncertainties. The fitted short glitch template is not 
removed from the data (see discussion on the impact later). 

4.1.2. Data flagging 

Data in the vicinity of the peak of the long or slow events for 
which the fitted long glitch template is above three sigmas of the 
noise are flagged and are not used for the processing (and for the 
next event fitting). For short glitches, data for which the short 
glitch template is above 0. 1 sigmas are flagged. 

At the end of the glitch detection and subtraction, a match 
filter optimized to detect 2 second exponential decays is applied 
to the data in order to find events which are left in by the method 
or events for which the fit failed. This introduce an extra flagging 
to a: 0.1% of the data. Additional glitch cleaning and use of 
different adaptive filters were found not to improve data further. 

4.1.3. Sky subtraction 

The sky signal is estimated at each iteration on each ring, de- 
tector by detector independently using the redundancies of mea- 
surements. Third-order spline coefficients are fit to nodes sep- 
arated by 1.5 arcmin to account for sub-pixel variations which 
could otherwise be detected as glitch signal. Flags determined 
at previous iterations are used to reject data, and estimated long 
glitch templates are removed prior to the sky signal estimation. 

Special sky estimation and glitch removal was required when 
scanning through or near bright objects, including planets and 
bright sources in the galactic plane to account for systematic 
errors in the subtraction of the sky signal that could be falsely 
flagged as glitches. These errors are due to two main effects; 
slow cross-scan variations of the pointing within the ring scan 
time, pointing uncertainties of the order of a few arcsec (Planck 
Collaboration VI 2013), and sub-pixel variations which are not 
captured by the spline coefficients. 

For rings with large signal in the galactic plane, we add to the 
noise variance the square of a term of signal reconstruction er- 
ror proportional to the signal (1% of the signal for submillimeter 
channels and 0.5% for cosmological channels). This automati- 
cally increases the threshold used for glitch detection at bright 
regions of the sky. Also, long glitch templates are not removed 
for high sky signal. The correction of the threshold is effective 
for all rings for submillimeter channels (545 and 857 GHz), and 
only a few percent of the data close to the galactic center for 
cosmological channels. Details of this treatment is described in 
detail in Planck HFI Core Team (2011). Furthermore, a special 
treatment is performed to process the data corresponding to the 
scanning of large amplitude planets. The motions of planets dur- 
ing one hour on the sky, which are mainly cross-scan, as well as 
the cross-scan pointing uncertainties are accounted. We fit the 
amplitude of the planet at each scan and rescale the averaged 
signal estimate on the ring with the estimated amplitude coeffi- 
cient to reconstruct the sky signal, and remove the best estimate 
of it from the input data. This neglects main lobe beam asym- 
metries. However, we find that it is a good approximation at this 
stage given the small cross-scan shift of the planet signal during 
a ring time. 
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Iterations between sky signal and glitch detection and sub- 
traction are performed up to the convergence to a sufficient ac- 
curacy ring by ring. For most of them, a total of six iterations are 
necessary. Details of the algorithm are described in Planck HFI 
Core Team (2011). 

4.1.4. Template estimation 

Long and short glitch templates are fit to stacked (normalized 
and averaged) large amplitude events of typically 1000 sigmas of 
the noise. The stacked data are produced using the median value 
for each bin corresponding to the sampling period index after 
the maximum. The median is used to avoid significant contribu- 
tion from other glitches. Figure 25 shows the difference between 
long glitch stacks and the best fit template with four exponen- 
tials in different energy windows for one bolometer (143-3a). 
The glitch template is determined using events with amplitudes 
between 1200 and 2400 noise sigmas. The difference between 
the stacked glitch data and the fitted template is small, but not 
compatible with zero reflecting the limit of the summed expo- 
nential model. The difference is systematic and of the order of 
1%, which is much smaller than the fitting errors described in 
Section 5.2.2 (see Figure 30 in particular). So we expect that the 
impact of this inaccuracy on the final results is negligible. The 
difference between the stacks computed over different ranges of 
glitch amplitudes, chosen between 300 and 5000 sigmas after 
correcting for the non-linearity of the relative glitch tail ampli- 
tude described in 2.2, and the glitch template is as small as for 
the window in which the fitted template is determined. Thus, we 
do not see evidence of significant variations of the shape of the 
template with the amplitude of the glitch over a wide range. The 
template used for the analysis has been estimated on a limited 
range of energies to avoid mismatch due to the non-linearity ef- 
fects. For small glitches, with amplitudes < 300 sigmas, there is 
a non-negligible difference between the stacked glitches and the 
fitted template. This is attributed to strong selection biases, e.g., 
Eddington bias, that affects the stacked signal, and not to real 
variations of the glitch template. This bias is clearly observed for 
stacked samples before the events which show a different mean 
value than a few tenth of second after the events. A systematic 
low frequency component is clearly visible on the stacks. 

4.2. Results 

The fraction of the data which are detected by the deglitching 
method as contaminated by glitches, and are not used for anal- 
ysis goes from 6% for detectors with the lower amound of long 
glitches, to 20%. Figure 8 of Planck HFI Core Team (2011) 
shows a segment of data for a bolometer with high rate of glitch 
before and after glitch template subtraction (this was using an 
older version of the method, but changes are too small to be ap- 
preciated in the figure). This illustrates the necessity of removing 
long glitch templates from the data for an accurate determination 
of the cosmological signal. 

The subtraction of long and slow glitch templates allows to 
reduce very significantly the contamination of cosmological and 
astrophysical signal. Figure 26 compares the power spectra of 
noise, orthogonal to the direction of projection on the sky, from 
raw three-point filtered data, with and without glitch template 
subtraction for all detectors from 100 to 353 GHz which are used 
for cosmology. Both data sets were flagged in the same way after 
glitch detection, and gaps are filled with white noise in addition 
to the result of gaussian smoothing of the data. Power spectra are 



computed on data chunks of size ss 100 rings which are averaged 
to cover the nominal mission. Estimation of the sky signal has 
been previously removed ring by ring for each bolometer using 
redundancies as described in Section 4.1, and after subtracting 
the averaged signal for each ring to avoid filtering the noise on 
time-scales larger than the ring scale. Before correction, power 
spectra are highly contaminated by glitches at frequencies be- 
tween 0.002 Hz and 2 Hz. The two-second tail is responsible of 
the excess between 0.01 and 0.1 Hz, while the intermediate time 
constant of 50 ms produces excess around 0.4 Hz, explaining 
why PSB-a have a higher relative excess at that frequency. We 
observe a very important reduction of residual noise after tem- 
plate subtraction. Some detectors have very few contamination 
from long glitches, this is the case for 143-5 and 143- la which 
correspond to the two lowest spectra before correction which are 
seen on the top right figure. The improvement after deglitching 
is effective but small for those detectors as expected, and their 
corrected power spectra are very close to the fundamental limit 
of the noise. The power spectra for those detectors are then in- 
dicative of the power spectra of noise without glitches. Power 
spectra for other detectors seem to converge to this limit after 
template subtraction. Accounting from the fact that the 1/f part 
of the spectra varies in amplitude from bolometer to bolometer, 
we can conclude that residual contamination from glitches is be- 
low the noise at all scales for all detectors. As mentioned in 3.8, 
one of the 100 GHz bolometer (and only this one) has an ex- 
tra family of glitches which is not accounted by our deglitching 
method. Some of these events are detected as "long" in the pro- 
cessing, and then the long glitch template is wrongly subtracted 
(those events are shorter than the long events). But the matched 
filter, designed to detect long glitch tails (either positive or neg- 
ative) allows to flag segments of data when this happens, lim- 
iting the effect. The power spectrum from this bolometer is the 
one with the largest excess at frequencies around 0.1 Hz. This is 
still a small effect which is accounted in the total noise budget 
(Planck Collaboration VI 2013). 

The comparison of power spectra of noise orthogonal to the 
sky projection is indicative of the quality of the correction from 
glitches but cannot be used to derive very accurately the contri- 
bution to the maps. The reason is that part of the noise is filtered 
due to the "noisy" estimation of templates by the fitting proce- 
dure which are then removed from the data. This is not the case 
for the signal projecting on the sky since the method iterated be- 
tween sky signal and glitch detection, template fitting and sub- 
traction. The contribution of glitches to the total power after data 
reduction at the map level or even at the ring level can be studied 
with simulations. This is done in the next section (5). 

As already shown, long glitches are in coincidence between 
bolometers from the same PSB pair. Consequently, glitches be- 
low the detection threshold or unsubtracted tails are expected to 
induce residual extra correlation between bolometer timestreams 
from the same pair. In order to evaluate the contribution from 
residual glitches we have computed the cross-power spectra be- 
tween bolometers. This is shown in Figure 27 between detec- 
tors at 100 GHz after averaging cross- and auto-power spectra 
computed on data segments of a: 100 rings over the full nomi- 
nal mission. Sky signal has been removed and data have been 
processed in the same way than previously described for auto- 
spectra. Bolometers from the same PSB pair have an extra corre- 
lation above 0.2 Hz which is of the order of 2 to 7% of the white 
noise level in the spectrum for 100 GHz bolometers. The corre- 
lation degree varies from 2 to 3% for PSB pairs at 143 GHz, 1 to 
4% at 217 GHz and 2 to 4% at 353 GHz. Simulations presented 
in the next section, show that residual noise excess is of the order 
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Fig. 25. Stacking of long glitches of bolometer 143-3a in different energy windows after normalizing each event to the max ampli- 
tude. The solid line correspond to the long glitch template obtained after fitting four exponentials using the data shown in the last 
window. The bottom plots show the difference between the stacked data and the template. We observe a good agreement between 
the template and the stacked data within typically 1%. 



of 5% (Figure 30) given our modeling of glitches (Section 3.2), 
and, combined with the fact that about half of glitches have the 
same amplitude in PSB-a and PSB-b (Section 2.7), we should 
have roughly 5% correlation in the noise spectrum. This is higher 
but comparable to the 3% observed for this pair. The difference 
could be due to uncertainties in our modeling and model ex- 
trapolations. Very importantly, we observe more correlations for 
pairs for which both detectors have less detected long glitches; 
e.g., the 100-4a/4b pair for which both detectors have the lowest 
glitch rates, as shown in Figure 9, and which have 7% cross- 
correlation (as seen in Figure 27), the highest observed com- 
pared to all pairs. This is entirely explained by the fact that both 
detectors have then more glitches below the detection threshold, 
since the counts for long glitches are similar between detectors 
but are essentially scaled in energy as described in Section 3.2. 
As opposed, we observe that pairs for which the faint-end break 
in counts is above the detection threshold have significantly less 
correlations above 0.1 Hz. 

We have clearly shown that glitches below the threshold are 
responsible from the extra-correlation of noise by averaging the 
two timestreams from two bolometers in the same pairs, and then 
by detecting events above 3.2 sigmas of the noise in the com- 
bined data. This allows to detect coincident glitches with am- 
plitude lower by a factor V2 for about half of the events. This 
allowed to divide the correlated noise by a factor of 2 with only 
0.1% of extra flagging of the data. This shows without ambiguity 



that extra-correlation between two bolometers in PSB pairs are 
due to undetected long glitches below the threshold. 



5. Impact of glitch residuals on final results 

5.1. Simulations 

In order to estimate residual contamination from glitches after 
processing, we have performed simulations of glitch timestream 
incorporating the important glitch properties we measure in data. 
Simulations include the following features: 

- Glitches are generated with a Poisson distribution with a sub- 
sample resolution. 

- Population spectra found in data for each population by com- 
bining bolometers, and the model explained in Section 3. In 
particular, for long glitches we use the measured counts for 
one of the submillimeter channels, which is very sensitive to 
long glitches, and for which we observe the break at low am- 
plitude with a very high significance. For this bolometer, we 
observe rising differential counts for events very close to the 
threshold. The counts are extrapolated at lower flux with a 
constant law (this has a negligible effect). Long glitch spectra 
are rescaled for each bolometer by the surface of the silicon 
die and energy is rescaled using a factor fitted on data. Short 
glitch populations are derived using one reference bolome- 
ter after correcting from the effect of non-linearity on the 
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Fig. 26. Power spectra of the component of noise which do not project on the sky, for each bolometer at frequencies between 
100 GHz (top left panel) to 353 GHz (bottom right panel). Black curves correspond to the power spectra without subtraction of 
templates for PSB-b and SWB, and the blue curves are for PSB-a. Power spectra are rescaled in amplitude such that they all match 
at 20 Hz. The red curves are with template subtraction. Power spectra are computed over « 100 rings and are averaged over the 
nominal mission. Sky signal is previously subtracted from data as described in the text. This creates the narrow dips in the spectra 
at the harmonics of the spin frequencies. Data have been previously three -point filtered and flagged after glitch detection. The large 
lines on the spectra, unsubtracted at this stage, are induced by the 4 Kelvin compressor and are described in Planck Collaboration 
VI (2013). The decrease of the spectra above 10 Hz is mainly due to the three point filter applied to data. Glitches contribute to the 
power spectra between 0.002 Hz and 2 Hz. Power spectra of residual contamination are significantly reduced after glitch template 
subtraction at frequencies below a few Hertz. There is much less dispersion between power spectra from different detectors after 
correction - the dispersion becoming less than a factor of two - indicating that remaining contamination from glitches is not a 
dominant effect. 



spectra by the approach discussed in Section 2.5, and low 
amplitude glitches below 20 sigmas of the noise are extrap- 
olated by a power law. The spectra are rescaled in energy to 
match the counts measured in data for each bolometer. Slow 
glitch population spectra are simulated using the measure- 
ment from one of the bolometers and extrapolated using a 
power law at low energy. Slow glitches are simulated only 
for PSB-a bolometers. As for other populations, the slow 
glitch spectra are rescaled in energy to match the measured 
spectra. 

- The non-linearity of the slow tail of long glitches relative to 
the peak amplitude as measured in data, see Section 2.2 for 
details. 

- The temporal shape of glitches is simulated using templates 
for the slow parts (after * 20 ms), correcting analytically the 
effect of the sampling average and of the three-point filter 
on the amplitude of the exponentials, as those are applied to 
data before estimation of templates. The fast part of glitches 



is simulated with a full resolution, and then averaged over 
the sample period using a simple boxcar average. 

Simulations do not include the apparent scatter of the slow tail 
amplitude of long glitches, apart from the intrinsic scatter due 
to the variations of the arrival time in the sample period reduc- 
ing the amplitude in a different way. Those fluctuations does not 
affect strongly the performances of the method since the slow 
part of glitches is fitted before subtraction, and so those am- 
plitude variations are absorbed by the fit. Moreover, we used 
a fixed template for short glitches, which is not realistic since 
we observed scatter in the slow tail, and we observed two sub- 
categories of short glitches. The effect of using a single, fixed 
template for short glitches should be very small since the short 
glitch template are not removed. We did not include the 2 sec- 
ond exponential decay in for short glitches in simulations since 
uncertainties in its amplitudes are important. Nevertheless we 
have evaluated that the contribution to the power spectrum in 
timestreams should be below 0. 1 % of the noise power spectrum 
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Fig. 27. Cross- and auto-power spectra between a few detectors 
at 100 GHz. Spectra are computed in the same way than de- 
scribed in Figure 26 caption. Black curves correspond to some of 
the auto-spectra shown in Figure 26. Data have been previously 
rescaled such that auto-spectra match at 20 Hz, and such that 
the white noise level approximately 1 . Blue curves correspond to 
cross-spectra between detectors from the same pair. Red curves 
are for detectors from different pairs. The extra correlation ob- 
served for PSB pair bolometers above 0. 1 Hz at the level of 2 to 
7% of the white noise spectrum is due to residual glitches below 
the threshold. Some residual correlation at the level of «0.1% 
at m 0.2 Hz is observed between other bolometer pairs which 
are spatially close to each other might result from imperfect sky 
signal subtraction due to pointing drift, as the error is coherent 
between different detectors nearby. The low frequency noise is 
highly correlated between detectors and is due to thermal fluctu- 
ations of the focal plane. 



given the population spectrum of short glitches. The effect of 
non-linearity of bolometers and of the ADC is not simulated. 
On data, the non-linearity has the effect of reducing the fast part 
of the glitches, which are short essentially since they dominate 
at very high amplitudes, relative to the slow part. We do not ex- 
pect this to degrade the performances of the method since the 
slow part of glitches since the template fitting starts a few sam- 
ples after the peak, in parts of the signal which is not affected by 
non-linearity. Finally, we did not implement the correlations of 
glitches between PSB pairs. The correlation of noise is not tested 
at this stage, and we treat each bolometer independently. 

To the glitch generated signal, we added gaussian noise con- 
taining a white component and a low frequency component de- 
scribed by a power low fitted to the data (index and /knee)- We 
also added a timestream of pure signal obtained by scanning a 
generated CMB map, as well as galactic dust and point source 
maps from the Planck Sky Model (Delabrouille et al. 2012), us- 
ing the pointing solution derived for the data. The constructed 
signal timestream is performed interpolating extracted signal 
from the map to limit sub-pixel effects (see Reinecke et al. 2006, 
for the methodology). We filtered the simulated data using the 
same three-point filter that for the data. The differential counts 
recovered from the simulations after processing is shown in 
Figure 28 for bolometer 143-2a. We observe a very good agree- 
ment with the spectra recovered in Planck bolometer data shown 
in Figure 1 1 . 



Fig. 28. Differential counts for the three populations of glitches 
for one of the detectors with high rate of glitches (143-2a) mea- 
sured after analysing simulated data in the same way than actual 
data. This can be compared to the results on real data shown in 
Figure 1 1 for the same bolometer. Dashed lines are shown for 
reference and are identical than the lines shown for real data. 
We find an almost perfect agreement between simulations and 
real data, except for very high amplitude glitches > 10 4 sigmas. 
Those glitches are affected by the bolometer non-linearity in real 
data, and we did not simulate the effects of non-linearity in sim- 
ulations, explaining the discrepancies in the counts. It also ap- 
pears that we simulated slightly less very high amplitude slow 
glitches than in data. This concerns only ~ 1 event per day. 

5.2. Error estimation 

5.2.1. Evaluation of signal bias by the deglitching procedure 

Due to the high signal to noise of HFI data, the sky signal sub- 
traction is a critical part of the deglitching procedure. Errors in 
the sky signal estimation could easily induce spurious detections 
of glitches and errors in the template subtraction which would 
then correlate to the signal. This could bias the sky signal estima- 
tion for two main reasons; first by flagging as bad data slightly 
more on average when the sky signal is higher (or lower); sec- 
ond, by subtracting slightly more glitch templates signals when 
the fluctuations of the sky signal is of a given sign. 

We have verified the absence of bias of the signal with the 
help of simulations. Indeed, we have computed signal rings by 
projecting pure signal timestreams used for simulations in one 
side, that we write rn;(f>), where p is the ring pixel and i is the 
ring number, and in another side we have projected the simulated 
observed data from which we have subtracted estimated glitches 
templates after applying the deglitching procedure as performed 
for real data, and also after rejecting data detected as contami- 
nated (flagging). We write this last quantity r ,(/;>) for data ring i. 
We then computed the average binned power spectrum of ro;(p) 
as: 

1 N 

Po(q) = 777T 7! 7! roiik) * r 0i (k)\ (2) 
""l i=l keD„ 

where f denotes the transpose conjugate, ro,(&) is the Fourier 
transform of ro,(p) and we summed over N rings (with N of 
the order of 10 000). We also computed the average cross power 
spectrum between ro,(p) an d h(p) as: 

1 N 

P ^ r ° i(k) * m ^ (3) 
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Fig. 29. Constraints on the signal bias after deglitching. Black 
curve: power spectrum (k 2 Po{k) in the text with k taken at the 
center of each bin q) of ring data after projecting pure simu- 
lated signal timestream for one of the detectors at 143 GHz (see 
below), containing CMB anisotropics, galactic dust and point 
source signal, and using the actual pointing; For the other dis- 
played points: we computed the cross-spectrum between 1. the 
estimated signal on rings after subtracting long glitch templates 
to the simulated data and after flagging, and 2. the input sig- 
nal only averaged on the same rings. Crosses and diamonds 
correspond to the difference between this cross-spectrum and 
the pure input signal power spectrum shown in black, for posi- 
tive and negative points, respectively. The quantity displayed is 
k 2 (P c (k)-Po(k)) as described in the text. We averaged the power 
spectra of ten thousand rings. Simulations are performed for one 
detector with a high rate of glitches at 143 GHz (143-2a). In the 
absence of bias of the signal due to the deglitching we expect that 
the cross-power spectrum described above is an unbiased esti- 
mate of the input signal power spectra, and then the difference 
should be compatible with zero. We do not detect significant bias 
in any of the 100 logarithmically spaced bins, and we place an 
upper limit of » 5 x 10~ 4 in each bin at the scales relevant for 
the CMB analysis. 

In absence of bias of the signal introduced by any of the two ef- 
fects previously described, P c (q) should be an unbiased estimate 
of P (q): 

P {q) =< P c (q) >, (4) 

where < . > is the ensemble average over infinite number of 
realisations of noise and glitches. 

In practice, the quantities Po(q) and P c (q) - Po(q) are evalu- 
ated for simulations of the 143-2a detector for N = 10 000 rings. 
We selected rings with no strong galactic signal in the galac- 
tic plane to avoid giving too much weight to the galaxy, keep- 
ing about 95% of the rings. Results are shown in Figure 29 for 
logarithmically spaced bins. We observe that P c (q) - Po(q) is 
compatible with zero at all scales. We measure a higher limit 
of the bias of w 5 x 10~ 4 in individual bins at all scales rele- 
vant for cosmology, i.e. for k corresponding to 1 < i < 2000. 
The absence of significant bias due to the deglitching procedure 
is verified at the ring level, but this conclusion can in principle 
be drawn at the map level and hence at the CMB power spec- 
trum level. Also, we do not expect to observe differences for 
the CMB measured with other bolometers since the presence of 
bias has to do with the errors made on sky signal reconstruction, 
and not on errors made of the glitch template fitting and sub- 
traction. Nevertheless, we have performed the same exercise for 



a SWB bolometer at 143 GHz (143-5). We do not find any bias 
and placed an even lower limit of 2 x 10 4 due to lower rate of 
glitches for this bolometer. 

The situation is different for strong galactic signal and for 
high frequency channels as the detection threshold is increased 
with the signal amplitude and, even, long glitch templates are not 
subtracted for very high sky signal, as described in Section 4.1. 
We expect this to bias slightly the estimation of the sky signal in 
the positive direction since less positive glitch signal is removed 
in regions of strong sky signal. This effect is described in detail 
in Planck HFI Core Team (2011). In particular we observed an 
effect of the order of 4x 10~ 4 of the signal amplitude at 545 GHz. 
The effect on the beam response estimation is studied in detail 
in Planck Collaboration VII (2013). 

Because the sky signal is estimated using data from individ- 
ual rings, this estimation is a noisy estimate with an rms which 
is roughly reduced by 7 as compared to the rms of the noise in 
each measurement. This induces some weak correlations in the 
glitch detection between each circle in the ring, after the sky sig- 
nal has been removed from the data. This effect does not bias 
the estimation of the signal as was shown before, but correlate 
slightly the remaining noise, after flagging and template subtrac- 
tion, between subsets of data measured by taking half rings, at 
a level of « 1%. We have studied this effect in detail in Planck 
Collaboration VI (2013) (see Table 2 and Figure 32 of this ref- 
erence in particular), and it has been taken into account for the 
noise prediction. 

5.2.2. Glitch residual contamination 

We have evaluated for some bolometers the level of contamina- 
tion by glitches left in data after the processing at the ring level 
by comparing the power spectrum of the input simulated noise 
projected on rings with the power spectrum of processed simu- 
lated data, after removing the input sky signal from the data. This 
is shown in Figure 30 for simulations of data of the same two 
detectors as in the previous section: 143-2a containing a high 
rate of long glitches; and 143-5 with a low glitch rate but con- 
sequently more glitches below the threshold (see Section 3.2), 
after averaging the spectra over 10000 rings. We also compare 
the two power spectra previously described with the power spec- 
tra obtained without removing estimated glitch templates to the 
data but using the same produced flags. 

For bolometer 143-2a, we can see that glitches dominate 
over the noise on large scales, even after flagging, for k < 2 x 
10~ 3 arcmin -1 . We observe an important gain in the noise power 
spectrum after removing long glitch templates. Nevertheless, 
glitch residuals still contribute to the noise power at the level 
of 30% below k m 2 x 10~ 3 arcmin -1 and below 10% above 
k x 6x 10~ 3 arcmin -1 . We attribute the excess at low frequency 
to errors in the subtraction of templates. Undetected glitches be- 
low the threshold contribute significantly at higher frequencies, 
at the level of 5% of the power for this bolometer, and causes 
the observed excess for 0.03 < k < 0.1 arcmin -1 . For bolometer 
143-5, the contamination by glitches before template subtraction 
below k < 10~ 3 arcmin -1 is between 15 and 20% and of the or- 
der of 8% above. After glitch template subtraction, the contam- 
ination reaches the level of 8% at almost all scales. This value 
corresponds to the level of contamination by undetected glitches 
below the threshold, which is higher than for the 143-2a bolome- 
ter in agreement with our modeling (see Section 3.2). 

This study shows that the remaining contamination from 
glitches is below the instrumental noise level even for channels 
with the highest glitch rates. The evaluation of the contamina- 
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Fig. 30. Estimation of power spectra of noise residual after averaging data on rings for the 143-2a bolometer (left) and the 143-5 
bolometer (right). Power spectra are averaged over 10000 rings and are computed after removing input simulated signal to the 
simulated data and using the estimated flags. The dot-dashed curve is without subtraction of long glitch tails, the solid curve is with 
subtraction as done for real data, the dashed curve is for the same simulations but without glitches. Residual glitches are significantly 
reduced after template subtraction. 



tion at the map level and on the CMB power spectrum has not 
been performed in this paper and is postponed to the release of 
the data of the full Planck mission. Nevertheless, the contribu- 
tion of glitches to the errors on cosmological results (Planck 
Collaboration XVI 2013; Planck Collaboration XV 2013) has 
been accounted since noise estimation is based on survey differ- 
ence maps (Planck Collaboration VI 2013). 

6. Conclusions 

In this paper we propose an interpretation of the glitch signal 
observed in Planck-HFl bolometer data, we present details of 
the glitch subtraction and flagging method, and we evaluate the 
impact of residual glitches in cosmological data after processing. 

We observed three families of glitches, short, long and slow, 
where each is described by a sum of 3 or 4 exponentials. They 
are distinguished from each other by the relative amplitudes and 
time constants of the 4 terms. The rate of each glitch type is an- 
ticorrelated with solar activity indicating that they are all due to 
energy deposition by galactic cosmic rays. Using the response 
of the on board particle monitor and cryogenic stages during so- 
lar flares, we find that the spacecraft and surrounding instrument 
mass shields the 0.1K plate and bolometers from all particles 
with energies at least < 39 MeV. Short glitches have a response 
similar to the optical response. These short glitches result from 
direct particle energy deposition in the grid for the lowest ener- 
gies or in the Ge thermistor for the highest energies. 

Long glitches dominate the total glitch count rate, but are 
usually smaller in amplitude than that short glitches. Long 
glitches show a fast part with ~ 10 ms decaying time constant, 
followed by the sum of slow exponentials from ~ 50 ms to 2 sec- 
onds time constants. We find that the flux of the long glitches is 
the same for all bolometers, when we use the area of the Si die. 
Coincident counts in the PSBs paired in a single bolometer mod- 
ule in the flight data, as well as laboratory tests of spare bolome- 
ters on the ground using alpha sources to expose the absorber 
grid and Si die separately, and tests using 23 MeV proton beam 
line, support the conclusion that long glitches are due to parti- 
cles depositing energy in the Si die. We find a minimum energy 
deposition in the count distribution of long glitches that corre- 
sponds to a proton with an energy at the minimum of the proton 



energy loss spectrum traversing the Si die normal to the surface. 
This minimum deposited energy is above the noise threshold for 
about half of the detectors. Thus we can detect individually al- 
most all of the glitches from the bolometer data. We suspect that 
the fast part is due to hot ballistic carriers propagating from the 
wafer to the thermistor. The slower time scale response is from 
the thermal relaxation of the Si wafer. We find evidence of long 
glitches caused by galactic cosmic ray helium but no evidence 
of galactic cosmic ray electrons. Low energy showers of sec- 
ondary or delta electrons were suspected to be the cause of the 
high rate of glitches in HFI. However, the high rate in PSBs is 
measured in coincidence between the two bolometers. This sig- 
nature cannot be caused by a delta electron shower. However, 
there is evidence of delta electrons in the coincidence of short 
glitches in PSB pairs. In this process, the delta electrons are pro- 
duced when a high energy particle hits the absorber grid of one 
bolometer and produces coincidence glitch at high probability 
and at lower energy, in the neighboring grid. 

Slow glitches are not understood at this time and have not 
been replicated yet in the laboratory. These glitches have nearly 
the same long time scale behavior as long glitches and thus are 
expected to be from energy deposition somewhere else in the 
module. 

We detect high energy secondary showers from very high 
energy events both as coincidence glitch detection among many 
bolometers and heating of the 0. IK plate. The rate ~ 5/hour and 
energy of these events is consistent with impact of galactic cos- 
mic rays at very high energies >1 GeV/nucleon. 

We detect a second type of high energy event with a different 
time response than that of the secondary showers, a characteris- 
tic energy ~ TeV and at a rate that is too high to be due to cosmic 
rays. We currently do not have an explanation for cause of these 
events. 

The new glitch analysis presented here is an improvement 
in the one previously reported. The method iterates between sky 
signal estimate ring by ring independently, and glitch detection 
above 3.2 sigmas of the noise and subtraction. The sky signal 
is carefully estimated before subtraction using spline interpola- 
tion to avoid spurious event detection. For a better separation 
of glitches into different categories, we employ for each de- 
tected event a joint fitting of the three glitch types using tem- 
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plates. Long and slow glitch templates are subtracted from the 
data. Without this subtraction, residual glitch tails after flagging 
would dominate over the noise in a large range of frequencies. 
The glitch removal method improves the noise performance of 
HFI, we found that the residual contribution from glitches after 
template subtraction is below the noise at all frequencies, and 
has maximum amplitude at around 0.1 Hz. Realistic simulations 
including the modeling of glitches support this result, and show 
for the two bolometer studied at 143 GHz that residual glitch 
contribute to less that 20% of the noise spectrum on scales along 
the scan larger than w 2 degrees. At smaller scales the resid- 
ual contribution is smaller than » 5% and is induced by long 
glitches below the threshold. This is supported by the observed 
correlations between the signals from bolometers in a PSB pair. 

The rapid part of all detected events is flagged. 
Consequently, the amount of data rejected for the analysis 
varies from 8% to 20% depending on the bolometer. 

With the help of the simulations, we have evaluated potential 
biases introduced by the method on the sky signal estimation. 
We did not find evidence of bias on cosmological signal and we 
estimated an upper limit of 10~ 4 on the spectrum. 
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